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ABSTRACT 

We present a method to estimate the lensing potential from massive galaxy clusters for given observational X-ray data. The concepts 
developed and applied in this work can easily be combined with other techniques to infer the lensing potential, e.g. weak gravitational 
lensing or galaxy kinematics, to obtain an overall best fit model for the lensing potential. After elaborating on the physical details and 
assumptions the method is based on, we explain how the numerical algorithm itself is implemented with a Richardson-Lucy algorithm 
as a central part. Our reconstruction method is tested on simulated galaxy clusters with a spherically symmetric NFW density profile 
filled with gas in hydrostatic equilibrium. We describe in detail how these simulated observational data sets are created and how they 
need to be fed into our algorithm. We test the robustness of the algorithm against small parameter changes and estimate the quality of 
the reconstructed lensing potentials. As it turns out we achieve a very high degree of accuracy in reconstructing the lensing potential. 
The statistical errors remain below 2.0 % whereas the systematical error does not exceed 1.0 %. 

Key words. (Cosmology:) dark matter. Galaxies: clusters: general. X-rays: galaxies: clusters. Gravitational lensing: strong. 
Gravitational lensing: weak 



1. Introduction 

The core structure of galaxy clusters provides important cosmo- 
logical information. Based on numerical simulations, we expect 
the dark-matter distribution to follow a universal profile with 
characteristic gradients and a scale radius ( [Navarro et al.|199"7) l. 
Outside relatively small central regions, cluster density profiles 
should not be strongly affected by baryonic physics because of 
the long cooling times in the intracluster plasma. Cold dark mat- 
ter is expected to clump on virtually arbitrarily small scales. The 
level of substructure in clusters thus potentially constrains the 
nature of the dark-matter particles ( |Boylan-Kolchin et al.|2009 
|Gaoetal.|20rT] ). 

The ratio between the scale and the virial radii of galaxy 
clusters, dubbed the concentration parameter, is frequently ob- 
served to be substantially different than theoretically expected. 
In particular in strongly gravitationally lensing clusters, concen- 
tration parameters significantly higher than those found in nu- 
merical simulations have been detected ( [Broadhurst et al.|2008[ 
|Coe et ar||2012[ Fig. 14). It is fundamentally important to find 
out whether this discrepancy reflects insufficient understanding 
at the level of our theory of cosmological structure formation, or 
whether it is a combination of baryonic physics, selection effects 
and measurement biases that gives rise to this observation. 

A multitude of precise observational data on galaxy clusters 
is or is becoming available: Weak and strong gravitational lens- 
ing constrain the distribution of the total matter density projected 
along the line-of-sight. X-ray emission and the thermal Sunyaev- 
Zel'dovich effect constrain the density, temperature and pressure 
of the intracluster gas. Galaxy kinematics constrain the gradient 
of the gravitational potential, albeit in a fairly entangled way. We 



should expect to find the strongest constraints on the core struc- 
ture of galaxy clusters by combining all available types of data 
in a common and consistent way. 



Without equilibrium or symmetry assumptions, only data 
from gravitational lensing can be interpreted, while the inter- 
pretation of gas physics and galaxy kinematics requires at least 
equilibrium assumptions. Given such assumptions, however, all 
types of galaxy cluster data can be theoretically modelled on the 
basis of the gravitational potential. In this first paper of an in- 
tended series, we are focussing on X-ray emission, devising an 
algorithm to convert the observed surface brightness profile into 
a projected potential that can then directly be combined with data 
from gravitational lensing. Studies in progress, to be reported on 
in due course, will concern the thermal Sunyaev-Zel'dovich ef- 
fect and galaxy kinematics. The ultimate goal of our studies is 
a non-parametric method combining strong and weak lensing, 
observations of thermal gas physics and galaxy kinematics into 
one consistent model for the projected cluster potential. 



This paper is structured as follows: In Sect. 2, we develop an 
algorithm (see |Lucy" 1974 1994 1 for reconstructing the projected 
cluster potential from the X-ray surface brightness. Numerical 
tests, described in Sect. 3, illustrate how this algorithm per- 
forms under reasonably realistic conditions. Although we adopt 
a spherically symmetric cluster potential for this test, spherical 
symmetry is not a necessary condition for our algoiithm to work. 
The influence of a possible deviation from spherical symmetry 
is exemplified in Sect. 3.3. The results and our conclusions are 
summarised in Sect. 4. 
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2. Recovering the projected gravitational potential 
from X-ray surface brightness 

2.1. Basic relations 

The central object of our study is the Newtonian gravitational 
potential O. Gravitational lensing measures the projection 
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geometrically weighted by a combination of the angular- 
diameter distances Dis.is from the observer to the lens, to the 
source, and from the lens to the source, respectively. Different 
lensing observables characterise derivatives of i/' of different 
order Time-delay measurements constrain differences in ijj 
along different lines-of-sight, deflection-angle differences be- 
tween components of a multiple image constrain differences in 
the gradient of ij/, elliptical distortions constrain the curvature 
matrix of i/r, and flexion will hopefully soon constrain its third- 
order derivatives. Combining all available lensing observables 
thus naturally leads to a reconstruction of the lensing potential 
ijj which can be more detailed where observables sensitive to 
higher-order derivatives can be measured. 

At least in or near hydrostatic equilibrium, the density and 
temperature of gas in the lensing potential well are also fully 
characterised by the Newtonian potential. We begin with the hy- 
drostatic equation 

VP = -pV<D (2) 

and assume that the gas satisfies the polytropic relation 
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with an effective adiabatic index y. Assuming spherical symme- 
try, equation (|2]l is immediately integrated to give 
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where quantities with a subscript refer to an arbitrary radius ro 
which could, for example, be set to zero. For practical reasons, 
we introduce a cutoff radius r^^ > ro and fix the potential such 
that Ocut - <l>(rcut) = 0. 

For y > 1 and a density profile that decreases monotoni- 
callMM <l)(r) can be arranged to be negative for r < ro by the 
structure of Eq. (Hli. Therefore p remains positive and semi- 
definite. The quantity 

Po 
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appearing in Eq. Q is the squared sound speed at the cutoff ra- 
dius. For convenience, we introduce the dimension-less potential 



,^ = 



r-1 
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and obtain the gas density 



P=Po^ 
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The temperature of an ideal gas in thermal equilibrium with the 
potential ip is 
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' It can be seen below that both cases are fulfilled in our considera- 
tion. 



where m is the mean gas-particle mass and kB is Boltzmann's 
constant. 

Since the frequency-integrated emissivity due to 
bremsstrahlung is given by 



jx-Cp'T'^' 
it can be related to the potential by 
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For realistic effective adiabatic indices 1.1 < y < 1.2 
( [Finoguenov et al.]|2001 1, the exponent rj is quite a large num- 
ber, 10 < 7/ < 20. 

Equation (10 1, together with the fact that ordinary lensing 



effects are determined by second-order derivatives of the pro- 
jected Newtonian potential, suggests the following algorithm for 
combining X-ray and lensing data: 

1 . By deprojection of an X-ray surface brightness map Sx, find 
an estimate jx for the X-ray emissivity jx- How this could 
be done e.g. by means of Richardson-Lucy deprojection will 
be discussed below. 



2. Use Eq. ( 10 1 to infer an estimate 
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for the three-dimensional, scaled Newtonian potential. 
3. Project ip along the line-of-sight to obtain an estimate iff for 
the two-dimensional potential, which is proportional to the 
lensing potential and can thus directly be combined with es- 
timates of Ip derived from lensing. 

Since rj is large, l/rj <& 1, which is a most welcome property 
of Eq. ( fTT) : Inevitable fluctuations in the deprojected estimate 
jx will be considerably smoothed that way. 

2.2. Deprojection 

Different algorithms exist for the deprojection of two- into three- 
dimensional distributions. Without symmetry assumptions, such 
algorithms cannot be unique. Assuming spherical symmetry for 
simplicity, a three-dimensional function f(r) is related to its two- 
dimensional projection g(s) by 



g(s)^j dz/(V.92+z2) 



(12) 



where s is the projected radius and the coordinate system is cho- 
sen such that the z-axis points along the line-of-sight. At fixed s, 
we have zdz = rdr, allowing us to transform Eq. ( 12 1 to 



g(s) = 2 
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which will be more convenient for our purposes. We rewrite the 
last equation 



-g{s)= drK{s\r)f(r), K(s\r) 
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where the factor l/n was introduced to ensure that the kernel K 
is normahsed. 
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dsK(s\r) = 1 
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with respect to integration over s. 

Richardson-Lucy deprojection begins with the generaHsed 
convolution relation 



g(y) 
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where the integral kernel K relates the variables x and y. By 
B ayes' theorem, the inverse problem 



(17) 
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fix) = fdyK'(x\y)g(y) 
has the deconvolution kernel 

K'(x\y) = -^Kiylx) , 

provided the deprojection kernel K is normalised as in Eq. (fT4| 
and the functions / and g are normalised with respect to the inte- 
grals over their domains. rlThese normalisations are a necessary 
condition for the algorithm to converge. 

Since f{x) is unknown, so is the deconvolution kernel. 
However, given an estimate fi(x) for the function f(x), a cor- 
responding estimate gi(y) of the projection is 
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allowing us to estimate the deprojection kernel K' by 

fi(x) 
K'(x\y) = ^:^K(y\x) . 

giiy) 
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Inserting this expression into Eq. ( [T7| ) gives the updated estimate 
fi+\ix) 



^ f,(x) j dy ^K(y\x) 
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for the function f{x) given the convolved function g(y), which 
will later represent the data. Beginning with a reasonable guess 
for foix), the iteration Eq. (21 1 usually quickly converges. 

A regularisation term needs to be included in presence of 
noise to prevent overfitting. |Lucy ( 1994 1 showed that, provided 
giy) is normalised, the deprojection algorithm described by the 
iteration ( |2T| can be cast into the form 
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containing the functional derivative of 
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with respect to fi{x). Where //[/] is equivalent to a likelihood 
function, which is maximised in order to obtain the best possible 
solution. He suggested to augment i/[/] by the entropic term 
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ff(x)dx = 1, fg(y)dy = 1 



containing aprior;if(x), to suppress small scale fluctuations. The 
functional H is then replaced by 



H[f]^Q[f]^H[f]+aS[f] 



(25) 



with a parameter a controlling the influence of the entropic term. 
Since 



^^=-ln^-l, 
Sfiix) X 

the entropic term changes the iteration prescription to 
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with 
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provided f is also normalised. 

This procedure completes our algorithm: We specialise the 
general deprojection kernel Kiy\x) to the projection kernel Kis\r) 
defined in Eq. ( 14 1, identify x with the radius r and y with the 



projected radius s. The function gis) is replaced by the measured 
X-ray surface brightness profile Sxis)- Then, Eqs. (22 1, (27 1 and 
(28 I allow us to iteratively reconstruct an estimate jxir) for the 
three-dimensional emissivity profile jxir), including an entropic 
regularisation term comparing the estimate jxir) to a prior ;\f(r). 
Recall that Sxis) and jxir) are assumed to be normalised as well 
as the kernel Kis\r). The complete iteration including the en- 
tropic regularisation term reads 



^jx,i - jx,i 
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The deprojection begins with a first guess jx.oir) for the X- 
ray emissivity profile, the prior xir) against which the deprojec- 
tion is to be regularised, and a parameter a controlling the degree 
of regularisation. 

The choice of a constant prior, x, leads to a statistical bias in 
the estimates of the deprojected functions such that they appear 
flatter than they should ( Narayan & Nityananda|1986 Lahav & 
Gull|1989[[Lucy|1994| ). 

This issue can be addressed by selecting as a default solution 
a smoothed version of the obtained result. This approximation, 
known as floating default ( Horne|1985[ Lucy|1994 i, is built by 
adopting the following definition: 



X(r) 
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where Pir\r') is a normalised, sharply peaked, symmetric 
function of r - r' and fir') corresponds to jxir)- In our investi- 
gations, we decided to choose a (properly normalised) Gaussian 
form with smoothing scale L: 



ir-r'f 



Pir\r') ~ exp - 



L2 



For practical considerations, we are working with discretised 
data sets and therefore the above integral formulation has to be 
approximated by sums. 



Sara Konrad et al.: Joint reconstruction of galaxy clusters from gravitational lensing and thermal gas 



3. Numerical tests 

3. 1. Simulating X-ray observations 




Fig. 1. Simulated image of photon counts from a galaxy cluster 
with amass of 5 x 10'"^ h^^Mg, a redshift of 0.2 and an exposure 
time of 1000 s. The detected photons have energies in the range 
of 0.5 - 8 keV. The resolution of this image was lowered by a 
factor of 6 compared to the resolution of the simulated Chandra 
image for visibility. 



The goal of this paper is to show that the algorithm defined in 
the previous section allows to recover the projected gravitational 
potential of galaxy clusters from their X-ray surface brightness 
profile, assuming spherical symmetry and hydrostatic equilib- 
rium. For testing this method we simulate galaxy clusters for 
which we choose a flat standard ACDM model with Q„, = 0.3, 
Q.b - 0.04 and Q./^ - 0.7. For the dark matter, building up the 
cluster potential well, we use an NFW (Navarro et al. 1997| l den- 
sity profile 

P('-) = ._,S\ ,,,,,.. > (33) 



(r/r,)(l+r/r,)2' 

with the scale radius r, and the characteristic density of the 
halo. We choose the concentration parameter c - 5.0. The gas- 
mass fraction is set equal to the universal baryon mass fraction 
fb = Q.hl0.m- The gas is assumed to consist of 75% hydrogen and 
25% helium, both completely ionised, with an eff'ective adiabatic 
index of y = 1 .2. 

We approximated the virial radius ryir of the simulated clus- 
ters with ^200, which is defined to be the radius within which the 
cluster's mean density is equal to 200 times the universal critical 
density at the given redshift. 
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The gas density and temperature profiles are then calculated, us- 
ing Eqs. (J6]l and (J8]l. To obtain a temperature profile which drops 
to zero at a large radius, we choose a large cut-off radius for 



the gravitational potential of rcut = 100 ■ raoo- The frequency- 
dependent emissivity due to bremsstrahlung is given by 



;x(v) = Cp^r-"^exp - 



hv 



(35) 



The expectation value for the number of photons, emitted in 
a detectable energy interval [£0^ ^i] pei" unit volume and time 
then reads 

^^ — d(MM^^ (36) 






hv 



with the cluster's redshift Zci- 

In order to obtain an image comparable to observations we 
simulate the CCD as foUows: 

- We neglect the convolution of the image with the telescope 
beam: Each pixel is mapped to a unique solid-angle element. 
The physical area 6A imaged by one pixel is 



5A = se'Dl,^, 



(37) 



where Dang is the angular diameter distance of the cluster and 
66 the angular side length of one pixel, assuming them to be 
perfectly quadratic. 

We did not include any cooling effects on the intracluster 
plasma even though cooling may steepen the X-ray surface 
brightness profile near the cluster core. 
Any absorption of X-ray photons between the cluster and the 
telescope is neglected. 

The detector has a perfect quantum efficiency within a sharp 
energy interval. Given the photon counts of the cluster given 
by Eq. ( 36 1, a pixel centred on the radial coordinate s is ex- 
pected to collect 
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photons per second. Dium is the luminosity distance to the 
cluster and Agfj is the effective telescope or detector area. 
Since Eq. ( 38 1 inherits a conversion from photon energy to 
photon counts, only one factor of (1 -H Zd) appears. 
The limited energy resolution of the telescope is mimicked 
by choosing appropriate energy intervals in Eq. (36 1. 



For the exact properties of the CCD, we adopt the same charac- 
teristics as the Chandra Advanced CCD Imaging Spectrometer 
(ACIS). The detection energy range is set to 0.5 - 8keV. We 
choose (Ziup - /iiow) /^res = 15 energy intervals, since our 
method is sensitive only to total numbers of photons and no 
spectral information is needed. We calculate photon numbers to 
pixels by drawing Poisson deviates with the appropriate expec- 
tation value 6N for all 15 energy intervals. The mean energy 
(E\ - Eo)/2 is assigned to each photon and the sum of energies 
allotted to the corresponding pixel. We include statistical noise 
by adding a constant background such that approximately 15 % 
of the detected photons is due to background. 

The pixel width is taken to be 0.5 arcseconds and the expo- 
sure time is set to 1000 s. In this way we obtain an X-ray sur- 
face brightness which is then azimuthally averaged around the 
centre of the cluster and binned. This profile is used as an esti- 
mate for the X-ray surface brightness profile and supplied to the 
Richardson-Lucy deprojection algorithm described above. 

Figures [T] and |2] show the photon counts and the normalised 
surface brightness of the synthetic observation and the nor- 
malised surface brightness profile of a simulated Chandra im- 
age for one realization of a galaxy cluster with a mass of 
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5 X 10 /j^'Mq and a redshift of 0.2. With these characteris- 
tics and the cosmology given above, this cluster has a scale ra- 
dius of approximately r^ - 0.25 /i 'Mpc and a virial radius of 

Tvir = 1.2/i"'MpC. 



1 0.01 - 

CO 

■^ 0.001 h 



Surface brightness 



0.2 0.4 0.6 O.E 

s/Mpc/h 



Fig. 2. Azimuthally averaged and normalised surface bright- 
ness profile of a simulated galaxy cluster with a mass of 5 x 
10'"* /j^'Mq and a redshift of 0.2. The corresponding synthetic 
observational data is shown in Fig. [T] 

The reconstructed and normalised gravitational potential 
is shown in Fig. [3h together with the true potential, the clusters 
were modeled with. Despite the statistical fluctuations of the sur- 
face brightness profile supplied to the algorithm, the contribution 
of the background noise exceeds the real surface brightness pro- 
file at large radii (i.e. s > 0.8 /i^'Mpc) which then leads to an 
overestimation of the gravitational potential. This effect and the 
normalisation condition of the deprojection algorithm leads to 
a slight underestimation at smaller radii (i.e. s < 0.3 /i^'Mpc). 
However, as we are only interested in the lensing potential, ma- 
jor fluctuations are averaged out as seen in Fig. [3b. 

3.2. Testing the algorithm 

Next, we test the sensitivity of the potential reconstruction 
against changing certain parameters of the deprojection algo- 
rithm. We compare reconstructions of the X-ray emissivity ob- 
tained by assigning diff'erent weights a from [0.0, 0.4, 0.7, 0.9] 
to the entropic regularisation term, with a = 0.0 corresponding 
to no regularisation, and constant L = 0.3 h^^Mpc. 

However, strong regularisation causes the reconstruction to 
overestimate the signal because of broad averaging by means of 
the flattening default kernel, because the regularisation penalises 
curvature. This can be seen at large radii. 

A similar conclusion can be drawn from Fig.Hb, which con- 
tains the deprojected profile for different choices of the smooth- 
ing scale L from [0.1,0.3,0.7,0.9] /j-'Mpc while a = 0.4 is 
kept fixed. Due to the smoothing process performed by the float- 
ing default regularisation, a large fraction of the noise pattern, 
especially for large radii, is averaged out, which improves the 
convergence towards the expected result. 

As a final test, we calculate the uncertainty of our recon- 
structed lensing potentials. We obtain this uncertainty by boot- 
strapping, sampling the synthetic cluster data A^ = 200 times and 



applying our reconstruction algorithm. For each result of the re- 
construction the mean squared deviation from the true potential 
is calculated and then averaged over the number of bootstraps, 
giving the rms deviation: 



rms(s; i/r) - 



\N 



1 ^ (^"°'-'"(5) - cr 

n=l Vt 



(s)r 



'C'ri^^r 



(39) 



where quantities marked with a superscript 'norm' are nor- 
malised to reach zero at the maximum projected radius. These 
rms are shown in the first panel of Fig.lS] Since this rms incorpo- 
rates statistical as well as systematic errors, we obtain a relative 
deviation from the true lensing potential of about 2 % for large 
values of s. The relative rms increases with increasing radius due 
to the poorer signal. 

Since we also know the exact surface brightness of the clus- 
ter simulation, we can estimate the relative systematic error of 
the algorithm itself. Doing so, we bin the real surface brightness 
profile in the same way as we binned our statistical CCD images 
and obtain a reconstruction that does not inherit any statistical 
fluctuations or noise. We call this the best possible reconstruc- 
tion (/^ideai- Its rms with respect to the true lensing potential (/rtrue 
provides an estimate for the systematic error of our reconstruc- 
tion algorithm (shown in the lower panel of Fig. |5]l. This mean 
squared deviation increases slightly with the cluster radius, but 
always remains below 1.0 %. 




Fig. 5. From top to bottom: rms deviation of iff from i^tme, rms 
deviation of lAideai from i^true and rms deviation of i/'subhaio from 
tptma, according to Eq. ( 39 1. Calculated from 200 realisations of a 
modeledgalaxy cluster with a mass of 5x10'"^ hT^M^, asubhalo 
with mass 1 x 10'"^ h^^Mg respectively, and a redshift of 0.2. The 
blue lines represent a 2.0% and 1 .0% threshold. 



3.3. Deviations from spherical symmetry 

As a last test we show how deviations from spherical sym- 
metry affect our reconstruction algorithm. For this purpose we 
create a CCD image of a single subhalo with a mass of 1 x 
10'"^ H^^Mq at redshift Zd - 0.2, having a virial radius of 
Asubhalo ~ 0.7 /i^'Mpc. We combine this CCD image with the 
simulated image of the cluster shown previously, with a pro- 
jected distance of 0.5 /z^'Mpc between the two cluster centres. In 
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Fig. 3. (a) Reconstructed and normalised gravitational potential of a simulated galaxy cluster with a mass of 5 x 10'"* h^^M^ and a 
redshift of 0.2. The potential was reconstructed assuming a = 0.4 and L - 0.3 /i 'Mpc. (b) Corresponding lensing potential for this 
simulated cluster. 
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Fig. 4. (a) Comparison between different choices for the weights a of the penalty function in the reconstruction of the X-ray 
emissivity with constant smoothing scale L = 0.3 h~^Mpc. The expected result based on Eq. ( 10 1 is plotted for reference. Whereas 
a = 0.0 correspond to no regularisation. (b) Comparison between different choices of the smoothing scale of the regularisation 
function in the reconstruction of the X-ray emissivity and a = 0.4. The expected result is plotted for reference. 



a first order approximation we accumulate the surface brightness 
of these two cluster realisations, before we apply our radial aver- 
aging scheme, again assuming spherical symmetry. Even though 
the subhalo is clearly visible in the lower panel of Fig. l5] (black 
arrow), we find that the maximum rms deviation of the lensing 
potential does not exceed the 2.0% deviation achieved without a 
subhalo. However, the reconstruction is rather insensitive to this 
kind of perturbation because of the radial averaging being ap- 
plied to the input data. Further simulations show that the closer 
the subhalo lies to the center of the hosting galaxy cluster the 
less is the influence on the reconstruction, since the degree of 
asymmetry decreases. 



4. Conclusions 

This work is motivated by existing and upcoming observational 
data on galaxy clusters. A set of accurate tools to constrain fun- 
damental cluster quantities is available: weak and strong gravi- 
tational lensing observations constrain the line-of-sight projec- 
tion of the cluster potential. X-ray observations constrain the 
density and the temperature of the intracluster gas, while the 
Sunyaev-Zel'dovich effect constrains the gas pressure. A further 
constraint on the gradient of the gravitational potential can be 
found by means of galaxy kinematics. Our final goal is to com- 
bine all non-lensing information on clusters with existing lens- 
ing reconstructio n methods (e.g.(B artelma nn et al.|1996[ Bradac 
let al.|2005|[2006l|Cacciato et al.|2006nMerten et al.|2009| l. 
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In this first paper of an intended series, we have outHned a Narayan, R. & Nityananda, R. 1986, ARA&A, 24, 127 
non-parametric reconstruction method for the projected gravi- NavaiTo, J. R, Frenk, C. S., & White, S. D. M. 1997, ApJ, 490, 493 
tational potential of galaxy clusters from thermal X-ray emis- 
sion, assuming hydrostatic equilibrium and spherical symmetry. 
It is more cumbersome, but quite straightforward to extend this 
method to triaxial halos. This extension is now being addressed. 

Our algorithm was laid out in Sect. 2. Assuming that 
the cluster is near or in hydrostatic equilibrium, we de- 
rived an analytic relation between the frequency-integrated 
bremsstrahlung emissivity and the three-dimensional Newtonian 
potential, Eq. ([TOji. The algorithm deprojects the observed X- 
ray surface-brightness profile by means of the Richardson-Lucy 
method, converts it to the potential and projects that. 

Numerical tests show how this algorithm performs under 
reasonably realistic conditions and how sensitive it is to its pa- 
rameters. We used one representative cluster realisation to ob- 
tain a realistic simulation of observational data. A comparison 
of different values of the weight a and the smoothing scale L 
of the regularisation term suggested that there is no significant 
variation in the final reconstruction of the X-ray emissivity if the 
floating default regularisation function is adopted with reason- 
able values. 

Even though our simulated galaxy clusters have a rather 
smooth surface brightness profile, this technique can be applied 
without restrictions to less well-behaved observational data, 
e.g. strongly peaked emission in the cluster centre due to cool- 
ing effects. In such cases, the peaked centre could be masked 
and then passed to the Richardson-Lucy algorithm. The results 
would still be reliable due to the local character of the recon- 
struction scheme. 

Furthermore, we determined the systematic error of our al- 
gorithm by applying it to ideally smooth rather than discretely 
sampled data to be at most 1.0%. We finally estimated the 
combination of systematic and statistical errors of our cluster- 
reconstruction algorithm to be at most 2.0% for s ~ ryir. Both 
rms were shown in Fig. E\ 

We also addressed the problem of a galaxy cluster with 
an asymmetric surface brightness, e.g. containing a subhalo, in 
Sect. 3.3. Due to azimuthal averaging, our reconstruction algo- 
rithm turns out to be rather insensitive to such kinds of perturba- 
tions. 
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